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Abstract. Propagation of pulses in myelinated fibers may be described by appropriate solutions 
of spatially discrete FitzHugh-Nagumo systems. In these systems, propagation failure may occur if 
either the coupling between nodes is not strong enough or the recovery is too fast. We give an 
asymptotic construction of pulses for spatially discrete FitzHugh-Nagumo systems which agrees well 
with numerical simulations and discuss evolution of initial data into pulses and pulse generation at 
a boundary. Formulas for the speed and length of pulses are also obtained. 
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1. Introduction. Effects of spatial discreteness are important in many physi- 
cal and biological systems comprising interacting smaller components such as atoms, 
quantum wells, cells, etc. Examples are the motion of dislocations |Dj, and crystal 
growth and interface motion in crystalline materials J6j, the motion of domain walls 
in semiconductor superlattices Q, pL sliding of charge density waves , pulse propa- 
gation through myelinated nerves |2j], and so on. The mathematical study of spatially 
discrete models is challenging because of special and poorly understood phenomena 
occurring in them that arc absent if the continuum limit of these models is taken. 
Paramount among these phenomena is the pinning or propagation failure of wave 
fronts in spatially discrete equations. Physically, the pinning of wave fronts is related 
to the existence of Peierls stresses in continuum mechanics (l7) , relocation of electric 
field domains |Q and self-sustained oscillations of the current in semiconductor super- 
lattices [|, electric current due to sliding of charge density waves 0, saltatory 
propagation of impulses in myelinated fibers and its failure Q , etc. 

Mathematical understanding of propagation failure of wave fronts in spatially 
discrete equations experienced significant progress after a paper by Keener po| . In 
Ref. pp| , Keener used comparison principles to characterize pinning of wave fronts 
and their motion for spatially discrete reaction diffusion equations of the form: 

(1-1) Un,t = d{u n+ i - 2u n + Un-i) + f(u n ), 

where / is a bistable source term and d measures the strength of the coupling. Mod- 

) include the spatially discrete Nagumo equation for nerve 
conduction |l2j, M or crystal growth || , and the Frenkel-Kontorova model for motion 
of dislocations ||13| . More recently, a number of results on the existence of traveling 
wave fronts u n (t) = v(n — ct) with smooth profiles v have been established [l(J . 

These papers do not characterize precisely the propagation failure of a wave front for 
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els described by Eq. (1.1 



critical values of the control parameters. For a piecewise linear source function, an 
explicit description has been given by Fath using extensively the properties of special 
functions [ jllj . In the general case, we have found that propagation failure can be 
explained as a loss of continuity of the wave front profile in critical values of the con- 
trol parameter 0. Furthermore, the smoothness of the wave front profile just before 
the propagation failure occurs can be exploited to obtain an analytical description 
of wave fronts and their speed near critical parameter values [Q, |^] . This theory has 
been extended to spatially discrete reaction-diffusion-convection equations describing 
dynamics of domain walls in semiconductor superlattices in Ref. g|. 

These advances in the mathematical understanding of propagation phenomena 
have occurred for spatially discrete scalar reaction-diffusion equations. Similar phe- 
nomena occur in models of calcium release at discrete sites |2^] . The latter consist of 
scalar reaction-diffusion equation with spatially inhomogeneous source terms that are 
close to f(u) times a series of delta functions centered at spatially periodic sites. Com- 
paratively, little progress has been made understanding wave propagation and fail- 
ure in spatially discrete systems. Anderson and Sleeman Q have extended Keener's 
techniques to discrete reaction-diffusion systems modelled by the FitzHugh-Nagumo 
(FHN) dynamics jl^, [2(| . Hastings and Chen have proved the existence of pulse 
traveling waves for a myelinated nerve model with a Morris-Lecar type dynamics. 
They also comment the difficulties of extending their results to the FHN system. An 
attempt to understand the mechanisms of propagation failure in the FHN system has 
been carried out by Booth and Erneux ||. They consider slow recovery and very 
special limiting (small) values of the parameters characterizing the bistable source 
and the spatial diffusivity in the FHN system. Furthermore, they also impose partic- 
ular boundary and initial conditions. With these restrictions, they could study how 
a specific disturbance localized in one cell propagated to neighboring ones until the 
resulting front failed to propagate. No construction of pulses or formulas for their 
velocity were given. 

In this paper, we construct asymptotically pulse solutions of the spatially discrete 
FHN system describing nerve conduction through myelinated fibers. We also discuss 
how the pulses may fail to propagate. Our ideas could be extended to spatially discrete 
systems whose cell dynamics contains widely separated time scales corresponding to 
fast excitation and slow recovery variables. Among these systems, let us cite models 
for bursting behavior in pancreatic j3 cells [^0) or the much more difficult case of 
front propagation in voltage biased semiconductor superlattices Q. In the latter, a 
separation of time scales exists but it is not obviously included as a small parameter 
in the equations. In our presentation, we have chosen the FHN dynamics for its 
simplicity. This model has been widely used to understand issues that are obscured 
by technical complications in more realistic models of nerve conduction. We consider 
the following system of dimensionless equations: 

d(u n+ i - 2u n + M„_i) + Au n (2 - u n )(u n - a) - v n , 
u n - Bv n , 

n = 0, ±1, . . .. Here u n and v n are the membrane potential and the recovery variable 
(which acts as an outward ion current) at the nth excitable membrane site (Ranvier 
node). The cubic source term is an ionic current and the discrete diffusive term 
is proportional to the difference in internodal currents through a given site. The 
constants A and B are selected so that the source terms in the FHN system are O(l) 
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Fig. 1.1. Nullclines for the space independent FHN model with different values of a. 



for u n and v n of order 1, that the only stationary uniform solution is u n — = v n 
and the FHN system has excitable dynamics (A = 1, B — 0.5 is a good choice, see 
Fig. [O]). The constant e > is the ratio between the characteristic time scales 
of both variables. We assume e <C 1, that is, fast excitation and slow recovery. A 



dimensional version of Eqs. ( |l.2| ) and (1.3) was derived in the Appendix of Ref. H 
from an equivalent-circuit model of myelinated nerves. For background on similar 
models, see @ ||, [H], ||. 

We shall study pulse propagation in the spatially discrete FHN system (1.2) and 



(1.3) by asymptotic methods. At first sight, such a task is hopeless: asymptotic 
methods require a degree of smoothness at appropriate time or length scales, and the 
spatial variable n in these systems is discrete. However, we can use the separation 
between time scales in the FHN system to show that a pulse is made out of two 
'sharp' wave fronts separating regions of slow spatial variation. Wave fronts are smooth 
solutions of the continuous variable z = n — ct/e, and perturbative arguments apply 
straightforwardly to them. Thus the theory of wave front propagation for spatially 
discrete scalar reaction-diffusion equations plays an important role in our construction 
of pulses. 

Let U\ < U2 < U3 denote the three zeros of the cubic nonlincarity f(u) in 
Eq. (1.1). Ui and U3 are stable solutions for d = 0. A wave front is a solution 
of Eq. (1.1) with a smooth profile u n (t) = u(n — ct) moving at a speed c, such that 
m(=Foo) = U\ and u(±oo) = U3. If /(it) is odd about u = U% and d is sufficiently small, 
a stationary solution of Eq. (1.1) exists and therefore no wave fronts can propagate 
JlO| |. As the source term departs from this symmetric form, front propagation is 
made easier. In Ref. [fjj, we selected d — 1 and f(u) — F — Ag(u), where g(u) is 
odd about its middle zero and F is an external force that quantifies departure from 
symmetry. Notice that we can obtain Eq. (1.1) with d = l/A and / = (F/A) — 
g(u) after rescaling time. We found that wave fronts propagate for \F\ > F C) where 
F c > depends on A and the specific g(u) we adopt. Equivalently, we could set 
f(u) = —u(u — a)(u — 2) and use a — 1 as a control parameter instead of F. After 
the "external force" a surpasses a critical value sufficiently far from the symmetry 
point a = 1, stati onar y fronts may cease to exist and propagating wave fronts may 
appear. See Fig. 1.2 (a). Notice that the limiting case considered by Booth and 
Erneux ("slow capture near a limit point"), d — 0(a 2 ), a — * 0+, corresponds to 
considering the parameter region in the lowest corner in this Figure. In this region, 
propagation failure is the normal situation. We could alternatively fix a and set 
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Fig. 1.2. faj Critical values a c i and a cr as functions of d. (b) Critical w c i (a, d) and w cr (a, 
for a = 0.5. (c) Critical w c i(a, d) and w cr (a, d) for d = 0.5. 



f(u) = — w — u(u — a) (it — 2), using w as the control parameter. Then there are critical 
values w c i(a, d) and w cr (a, d) such that wave fronts fail to propagate if w c i < w < w cr ; 
see Fig. |l.2| (b). How does the parameter a affects the critical values of w7 Assume 
that a < a c i(d) for a fixed value of d, so that wave fronts propagate for w = 0. To 
compensate this effect, we need a small critical value w c i{a, d) of the parameter w. As 
a is departs more and more from a c i(d), larger and larger critical values w c i(a,d) are 
needed to return to the situation of propagation failure. A similar situation occurs 
with w cr (a, d). Th us t he critical values of w increase (in absolute value) as \a — a c \ 
increases; see Fig. 1.2(c). The effect of increasing the diffusivity d is to shrink the 
parameter range where stationary fronts exist. In fact, as d — > oo (continuum limit), 
the width of the pinning interval is conjectured to decrease exponentially fast to zero 
for certain nonlincarities |6|, ^3). Propagation failure can be understood as a loss of 
continuity of the moving front as appropriate critical parameter values are approached. 
Increasing the discrete diffusivity and deforming the source term sufficiently far from 
odd symmetry about its middle zero, both facilitate propagation of wave fronts 
Reciprocally, weakening the coupling between cells and diminishing the "external 
force" a — 1 helps finding propagation failure. 

For the spatially discrete FHN system, the description of wave propagation is more 
complicated. This also happens for the spatially continuous FHN system (D d 2 u/dx 2 
instead of discrete diffusion). Depending on the initial condition, stable wave trains 
or pulses may be approached as time elapses [15], Q. Pulses cannot be obtained 
for e larger than a critical value. For discrete diffusion, we can construct pulses 
provided e is smaller than a critical value e c (a,d), a is outside a certain interval 
(corresponding to propagation failure in the scalar case) and the initial condition 
is chosen appropriately. Our construction combines the theory of front depinning 
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developed in JjJ with Keener's asymptotic ideas developed for the FHN model 
with spatially continuous diffusion. Our results agree very well with direct numerical 



solutions of (tLa) and (1.3) 



The key ideas of an asymptotic theory for (1.2) and (1.3) in the limit as e — > 
are simple. First of all, pulses consist of regions where u n {t) vary smoothly with n, 
separated by moving sharp interfaces (fronts). In the first type of region, we may set 
e = d = and obtain a description of slow recovery. The sharp interfaces are wave 
fronts with smooth profiles, u n (t) = u(z), v n (t) = v(z), with z = n — ct/e. Then v 
is constant at each side of a front and the excitation variable u obeys the spatially 
discrete Nagumo equation, whose fronts we can characterize £|. A stable pulse is 
obtained when the velocity of the leading front is equal to that of the trailing front 
JL9| . This condition fixes the pulse width. Its violation or propagation failure of any 
of the fronts bounding the pulse result in propagation failure thereof. Notice that 
we use an analytic expression for the wave front velocity of the Nagumo equation 
valid as a is near its critical value for propagation failure. For small values of d, 
the front propagation range is narrow and the formula for wave front velocity holds 
for all appropriate values of a; see Fig. |l.2| (a). For larger values of d, the interval 
where propagation occurs is wide and we can use our approximation only for a close 
to its critical values a c \ and a cr . Outside these parameter ranges, the velocity of the 
Nagumo wave fronts should be calculated numerically. 

The rest of the paper is organized as follows. In Section we recall certain 
needed results on wave front propagation and failure for the spatially discrete scalar 
reaction-diffusion (Nagumo) equation. Section || contains our main theoretical ideas 
in this paper, with the asymptotic construction of pulses for the discrete FHN system. 
These ideas and our results are tested by numerically solving the FHN system with 
appropriate boundary conditions in Section ^. Comments on propagation failure of 
pulses in the FHN system are made in Section |[ Section ^ briefly discusses how a 
pulse may be generated by applying a temporary stimulus at one end of a fiber with 
finitely many nodes. The last Section contains our Conclusions. 

2. The spatially discrete Nagumo equation. We consider the equation 
du 

(2.1) — - = d(u n+ i - 2u n + Un-i) + u n (2 - u n )(u n - a) - w, 
as 

for some constant w and denote h{u,w,a) = u{2 — u)(u — a) — w. As long as 
min/i(u,0, a) < w < max h(u, a, w) this is a 'cubic' source having three zeroes 
Ui(w,a), i — 1,2,3, U% < U2 < U3. Wave front solutions joining U\ and U3 (the 
two stable zeros) exist. A theory of pinning and propagation of fronts for this type of 
equations has been developed in |t], ||]. We sketch below its implications for (2.1). 



Assume first w = so that the assymmetry of the source is controlled by the 
parameter a. For d fixed, there are values a c i(d) and a cr (d) such that: 

• The fronts joining u = and u = 2 are stationary if a c i(d) < a < a cr (d). No 
front propagation is possible. 

• Outside this interval, there exist traveling wave fronts u n (s) = u(n — cs) 
joining and 2. For a > acrid), increasing fronts move to the right and 
decreasing fronts move to the left. For a < a c i{d), fronts move in the opposite 
way: decreasing fronts move to the right and increasing fronts move to the 
left. 

The values a c i{d) and a cr (d) can be approximately calculated as follows. In a large 
lattice, we decrease or increase a from 1 till we obtain an stationary solution u n (a) 
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whose linear stability problem has a zero eigenvalue. See Fig. 1.2. 

Now we fix a and vary w. The asymmetry of the source is controlled by a and w. 
For fixed d and a, critical values w c z(a, d) and w cr (a, d) are found such that: 

• The fronts joining Ui(w, a) and Usiw^a) are stationary if w c i(a,d) < w < 
w cr (a,d). 

• Outside this interval, there exist traveling wavefronts u n (s) = u(n—cs) joining 
U\(w,a) and 1/3(111, a). For w < w c i, these fronts move to the left if they 
increase from U\ to U3, and to the right if they decrease from U3 to U\. For 
w > w cr , fronts decreasing from U3 to U± move to the left, and increasing 
fronts move to the right. 

To calculate w c i and w cr , we start by fixing a and finding a value, w = Wo, at which 
stationary solutions exist for a large lattice. We now decrease or increase w from this 
value till we obtain a stationary solution u n (w) whose linear stability problem has a 
zero eigenvalue. See Fig. [E^. 

For w near any of its critical values, we can use the following formula to predict 
the speed of the fronts for |to| > \w c 



"c\- 



ton\ 1 a \ ■ 1 ,yja(3(w-w c ) 
(2.2) c(a, d, w) ~ sign(u; — w c ) . 

The parameters a and (3, given by a = 't'n, P = \ 6u n (w c ) + 2(2 + a)]<^ [Q, ||], 
are functions of a, d and the critical value of w. In these formulas, <p is a positive 
eigenfunction of the linear stability problem for u n (w c ) with ^ <^>„ = 1 and u n (w c ) is 



a stationary solution of (2.1) with w — w c ||. If w is not close to its critical values, 
the speed c(a, d, w) should be calculated numerically. 

A peculiarity of the Nagumo equation is the scenario for front propagation failure. 
As we approach the critical values for a, w, or any other appropriate pinning control 
parameter, the front profiles become less smooth and a number of steps appear. In 
the limit as the control parameter tends to its critical value, the transition regions 
between steps become infinitely steep, the front profile becomes discontinuous and its 
velocity vanishes ^ [|. 

3. Asymptotic construction of pulses. As we will discuss below, an appro- 
priate initial condition evolves towards a pulse. In particular, we need to fix the 
parameters d > 0, a < a c i(d) (the case a > a cr (d) follows by symmetry) and e smaller 
than a certain critical value, e c (a, d). This latter condition also holds for the spatially 
continuous FHN system which has two pulse solutions (one stable and one unstable) 
for e < e c . These solutions coalesce at e c and cease to exist for larger e ^J. A 
pulse consists of regions of smooth variation of u on the time scale t, separated by 
sharp interfaces in which u varies rapidly on the time scale T — t/e. In the regions 
where u varies smoothly, we can set e = d = 0, thereby obtaining the reduced problem, 

(3.1) tt„(2 - Un)(a - u n ) - v n = 0, 

dv n 

(3.2) ——=u n -Bv n . 

These regions are separated by sharp interfaces (moving fronts) at which u n varies 
rapidly as u n (i) — u(z), v n (t) — v(z), with z = n — ct/e. There, to leading order, 

du 

(3.3) -°-^ = d i u ( z + 2u ( z ) + u ( z - !)] + u ( z ) ( 2 - u ( z )) ( a - u ( z )) - v , 

(3.4) -of = 0. 

az 



Thus v is a constant equal to the value v n (t) at the last point in the region of smooth 



variation before the front. Equation (3J3) has a wave front solution as discussed in the 
previous Section. We can now discuss different regions in the asymptotic description 
of a pulse: 



1. The region of smooth variation of u in front of the pulse, described by (3.1) 

(v n ), so that 

Ui{v n ) - Bv n , 



and (3.2). In this region, u n = U±(v n ), so that 

~~dt 



and initial data evolve exponentially fast towards equilibrium, u n = v n = 0. 
The pulse leading edge. Let v(t) be the value of v n at the last point of 
the region in front of the pulse. Eventually, v — > 0. At the leading edge, 
u n (t) — u(n — ct/e) is a wave front moving towards the right with speed 
C = c(a,d,v)/e measured in points per unit time t. We have the boundary 
conditions u{— oo) = Us(v) and tt(oo) = U\{v) for the monotone decreasing 
profile u(z) which satisfies ( |3.3[ ). It is convenient to call c_(u) = c(a, d, v). 
Eventually, C ~ c_(0)/e, and u n decreases from u n — 2 to u n — across the 
leading edge of the pulse. 
Region between fronts: u n — U^{v n ) and 

dVn tt I \ D 
— — = U 3 (v n ) - Bv n . 

at 

There is a finite number of points in this region. On its far right, v n — v — > 0. 
As we move towards the left, v n increases until it reaches a certain value V(t) 
corresponding to that in the trailing wave front. 



4. Trailing wave front: v n (t) = v(z) = V, and u n (t) = u(z) obeys (3.3) with 
boundary conditions u(— oo) = U±(V) and u(oo) = U^{V). This front in- 
creases monotonically with z and it moves with speed C — c(a, d,V)/e mea- 
sured in points per unit time t. It is convenient to denote c+(V) = c(a, d, V). 
We shall indicate how to determine V below. Clearly, if the pulse is to move 
rigidly, we should have c + (V) = c_(0) after a sufficiently long transient pe- 
riod. 

5. Pulse tail. Again u n = Ui{v n ) and dv n /dt = U\(v n ) — Bv n . Sufficiently far 
to the left, v n = u n = 0. 

The number of points between wave fronts of the pulse is not arbitrary: it can be 
calculated following an argument due to Keener for the spatially continuous case . 
Let r be the delay between fronts, i.e. the time elapsed from the instant at which the 
leading front traverses the point n = N to the instant when the trailing front is at 
n = N . Clearly, 



(3.5) 



v ^ dv 



The number of points between fronts, l(t), can be calculated as 



(3.6) I = - / c-(v(t))dt 



e 



t-T 



On the other hand, the separation between fronts satisfies the equation 



The three equations ( |3.5|) , fl3.6| ) and (3.7) can be solved to o btai n the three unknowns 
r, Z and V(f). (The function v(t) is determined by solving (|3.2| ) with u„ = f/i(u n ) in 
the region to the left of the leading front). 

After a transient period, v(t) — > and 7(i) — > V (a constant value), so that we 
have the simpler expressions 



(3.8) 
(3.9) 



v 



dv 



, U 3 (v)-Bv 
dl _ c_(0) - c+(V) 



instead of (3.5) and (3.7), respectively. The number of points at the pulse top is now 



(3.10) 



c_(0)t c_(0) 



i7 3 (u) - Bv 



This equation yields V as a function of L Then (3.9) becomes an autonomous 
differential equation for I that has a stable constant solution at I = I* such that 
c_(0) = c+ (7(0): At I = r, the right hand side of Q has a slope -[C/ 3 (7) - 
BV]d + (V)/c-(0) < 0. 

Recapitulating, for appropriate initial conditions, leading and trailing fronts of a 
pulse evolve until I reaches its stable value at which c_(0) = c+(V(l*)) and ( 3.1 0| ) 
holds. To compute I*, we first determine V* = V(l*) by using c_(0) = c + (V(l*)). 
Then we calculate r = r* (which does not depend on e!) from ( |3.8| ) and I* = 
c_(0)r*/e. Our construction breaks down if the number of points between fronts 
falls below 1. This yields an upper bound for the critical value of e above which pulse 
propagation fails: e c ~ c_(0)r*. 

The asymptotic length of the pulse tail is obtained by first calculating the time 
needed for v n to go from to V(l*) to the left of the trailing front: T = dv/[Ux(v) — 
v]. The tail length is then L = c_(0)T/e. 

4. Numerically calculated pulses. We shall compare numerical solutions for 
different representative values of d with the approximate pulses provided by our theory. 
As initial data, we have adopted our approximate pulses. We have also used hump- 
like profiles with compact support for u n (0) and v n (0). It is important that v n (0) = 
at the leading edge of u n (0) and to its right, and that v n (0) = V rs w cr (d) at the 
trailing edge, where v n (0) reaches its maximum. Had we chosen v n (0) = for all 
n, the u n profile would split in two pulses traveling in opposite directions as time 
elapses; see Fig. 



4.1 



The region between leading and trailing fronts of the pulse 
acquires its asymptotic shape quite fast, but the pulse tail is usually rather long and 
evolves slowly towards its final form. 

• d = 1. In this case, a c i = 0.996. We choose a = 0.99 < a c i and analyze 
first front propagation for the rescaled Nagumo equation (2.1). The front 
propagation thresholds for these values of d and a are w c i = 0.0038 and 
w cr — 0.0095. Fig. [h^(a) shows the speeds of leading and trailing fronts 
as functions of w, as predicted by Eq. ( |2.2| ). For w = 0, the leading front 
should move at speed c_(0) = 0.0093. The relation c+(V) = c_(0) yields the 
asymptotic value V* = 0.0133 at the trailing front joining U x (V*) = -0.00665 
to U 3 (V*) = 1.9933. The time elapsed between fronts is r* = 0.00652, 

Then our upper bound for the critical value of 
Choosing a smaller value, e = 0.000005, we obtain a 
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Fig. 4.1. Splitting of an initial profile into two pulses propagating in opposite directions for 
d = 0.1, a = 0.5 and e = 0.006. 



pulse speed C = c_(0)/e = 1869 points per unit time and a pulse width of 
I* = Ct * ~ 13 poi nts. Our numerical solution of the full FitzHugh-Nagumo 
system ( |l.2| ) and (1.3) yields a pulse speed C = 2000 and a width of 13 
points for e = 0.000005. The trailing front joins -0.006647 and 1.993 with 
V* = -0.0133. See Fig. fO|(b). Note that the relative error in the predicted 
speed C is 0.0655. Obviously, rescaling the speed to C = c_(0)/e amplifies 
the error in our predictions. We have not been able to observe pulses for 
e > 0.0000076 which is smaller but not far from our estimation e c = 0.000064. 
Let us choose now a = 0.5 which is far from a c /. Then, w c i — 0.3194 and 
w cr = 0.3287. Eq. predicts c_(0) = 0.09983 whereas the trailing front 

joins —0.316 to 1.71 at V* — 0.6. If e = 0.01, the speed and width of the pulse 
are C — 9.983 and I* = 0.351C ~ 4, according to our theory. Numerically, we 
observe C = 64.7 and I* = 25. The source of these large errors is the value 
c_(0) = 0.09983 predicted with formula (2.2). If we replace this value by the 
numerical front speed calculated directly, c_(0) = 0.673, we obtain C = 67.3 
and I* ~ 24 points, which fit better the numerically observed values. 
d = 0.1. In this case, a c i = 0.567 and we shall choose a = 0.5 < a c i. Let us 
first analyze front propagation for the rescaled Nagumo equation ( |2.1|). F or 
these values of d and a, we obtain w c i = 0.0307 and w cr = 0.6175. Fig. |4.3j(a) 
shows the predicted speeds of leading and trailing fronts as functions of v, as 



given by formula (2.2). For w = 0, the leading front should move with speed 
c_(0) = 0.075662. At the trailing front, c+(V) = c_(0) yields V* = 0.648 and 
Ui(V*) = -0.33328 and U 3 (V*) = 1.666. The time elapsed between fronts 
is t* = 0.39266, which gives I* = 0.0297/e. Our bound for the critical value 
of e is e c = c_(0)t* = 0.029. Selecting e = 0.003, we predict C = 25.22 and 
I* ~ 10 points. Direct numerical calculations yield a pulse speed C = 26.38 
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Fig. 4.2. (a) Predicted speeds for wave fronts of the Nagumo equation {2.1) with d = 1, 
a = 0.99. The horizontal line marks the condition c+(u>) = c_(0), thereby graphically yielding 
1 1 = V* . (b) FHN pulse for e = 0.000005 . (c) Predicted and numerical speeds for wave fronts of Eq. 
(2.1) with d = 1, a = 0.5. The horizontal lines mark c_|_(«j) = c_(0). (d) FHN pulse for e = 0.01. 



and a pulse width about 10 poin ts. The trailing front joins —0.3269 to 1.675 
with V* = 0.6578. See Fig. fO](b). We have not been able to obtain pulses 
for e > 0.007, which is four times smaller than our upper bound of 0.029. 
• d = 0.01. In this case, a c i — 0.195 and we shall choose a = 0.1 < a c /. Let us 
analyze first front propagation for the rescaled Nagumo equation (|]l]). The 
front propagation thresholds for these values of d and a are w c i = 0.0136 and 
w cr = 1.0784. Fig. 4.4(a) shows the predicted speeds of leading and trailing 
fronts as functions of w according to Eq. (2.2). For w — 0, the leading front 
should move with speed c_(0) = 0.052. Then the trailing front has have 
V* = 1.092 corresponding to c+(V) = c_(0), and it joins Ui{V*) = -0.6 to 
UziV*) = 1.4. The time elapsed between fronts is r* = 0.748 and the pulse 
width, I* = 0.297/e. Our bound for the critical value of e is e c = c_(0)r* = 
0.058. Selecting e = 0.001, we predict C — 52 and I* = 39 points. Numerical 
observations yield C — 77.7 (a relative error of 0.3) and a pulse width of 59 
points. Furthermore, the trailing front joins —0.59 to 1.4 with V* = 1.095. 
See Fig. [0|(b). Again the observed errors in the pulse speed and width are 
due to errors in the prediction of c_(0) given by formula ( |2.2p . Replacing 
this value by the numerically computed front speed c + (0) = 0.078, we obtain 
C = 78 and I* ~ 58, better fit to the real values. 
Let us now describe the situation for other values of d. Our asymptotic theory 
agrees with the numerical results provided e is sufficiently small, but the velocity of 
the Nagumo wave fronts should be either approximated by Eq. (2^) or calculated 
numerically depending on how close to zero w c i happens to be. For d < 0.01, the 
length of the intervals where fronts of the Nagumo equation pro pagate is very small. 
Then the front speeds are always very small and given by Eq. (2.2) with great accuracy. 
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Fig. 4.3. (a) Predicted speeds for the Nagumo equation (2.1) with d = 0.1 and a = 0.5. The 



horizontal line graphically yields V* such that c_|_(V*) = c_(0). (b) Profiles of the FHN pulse for 
e = 0.003 (c) Trajectories of one point, uo(t), Vo(t), as the FHN pulse propagates through it. 



Our asymptotic description of the pulse agrees very well with numerical solutions of 
the FHN system. If d > 1, the spatially discrete FHN system can be approximated 
by its continuum limit. The length of the pinning intervals for the Nagumo equation 
is below 0.001, and the wave front velocities are essentially a correction of the wave 
front velocities for the spatially continuous Nagumo equation [^lj : 

(4.1) <=^{'- k 4 +o {i 

(A0 , 2U 2 (w )-U 1 (w)-U 3 (w) 

(4.2) c = 



(4.3) h = - 



V2 

fZ,e-&V{(s)Vr(s)ds 



124JZ o e-^(s)V^'(s)ds 
Here Vq is the appropriate wave front solution of the equation 
(4.4) Co 2 F " - V£ + V (2- V )(V -a)-w = 0. 

5. Propagation failure. Two facts may lead to propagation failure: a value of 
e that is too large or a <E (a c i(d), a cr (d)). 

Let us consider the first cause of propagation failure now. If e surpasses a certain 
critical value e c , recovery is too fast and a stable pulse cannot be sustained. This 
situation also occurs in spatially continuous FHN systems. In these systems, there 
exist two pulses (one pulse is stable, the other unstable) for e < e c , they coalesce at 
e c and cease to exist for larger e. In the discrete FHN system, the phenomenon of 
wave front propagation failure implies that pulses may propagate only if a < a c i{d) or 



a > a cr (d). As indicated by Eq. ( 3.1C ), the number of points between the two fronts 
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of the stable pulse decreases as e increases towards e c (a, d). Eventually the two fronts 
coalesce and it is not possible to propagate a stable pulse for e > e c (a, d). If we start 
with an appropriate pulse-like initial condition, we find the scenario of propagation 
failure depicted in Figs. 5.1 and [5^. For small d (d = 0.1), the variable v n ceases 
to be almost constant at the leading edge of the pulse and the distance between the 
two fronts diminishes. While v n ~ at the rightmost point of the leading front, 
v n ~ w > at the leftmost point. Thus, u n in this front decreases from U^{w) to zero 
as n increases. The value w increases with e and Uz{w) decreases. At the same time, 
the leading front speed diminishes as w increases till w surpasses the propagation 
threshold and the leading front stops. Since the back front goes on moving, the pulse 
vanishes, see Fig. 5.1. For large d (d = 1), a decremental pulse is formed. Its width 
and height decrease as it moves until it disappears, see Fig. |5.2| . Numerical simulations 
of the FHN system show that e c — * as a tends to either a c i(d) or a cr (d). 

Let us assume now that a G (a c i(d), a cr (d))- Then, the leading front cannot 
propagate with v n = v = 0. We need v n ~ v < w c i{a, d) < 0. However, in the region 
in front of the leading edge, v n and u n evolve towards 0, whereas we have u n > 
at the leading front. Thus dv n /dt = u n — Bv n > there and v n will increase until 
v n > 0, which contradicts our previous assumption. Thus we cannot have stable 
propagating pulses. Furthermore, there are no stationary pulses of the type we have 
discussed for this range of a: if v n = u n /B, the source u n (2 — u n ){a — u n ) — v n = 
u n (2—u n )(u n — a)—u n /B has only one zero, not three as in our construction. This does 
not preclude existence of other pulses, such as those corresponding to the homoclinic 
orbit in the phase space of the spatially continuous FHN system. However, we have 
not observed stable stationary pulses of this type in the spatially discrete FHN system. 



6. Pulse generation at a boundary. So far, we have considered the motion 
of a pulse (or its failure) in a sufficiently large myelinated nerve fiber. We have not 

12 




discussed how such a pulse might be created in a more realistic situation. Clearly, 
nerve fibers have finitely many Ranvier nodes and pulses are typically generated at 
the fiber boundary. Thus we are led to consider how a pulse might be generated by 
an excitation at a boundary and how the pulse propagates or fails to in a finite fiber. 
This problem was tackled by Booth and Erneux [^| using parameter values for which 
the FHN pulse fails to propagate. We shall now discuss different parameter ranges. 

Nerve fibers may have either a few Ranvier nodes (e.g., 20 for the central ner- 
vous system pq ]) or several hundred nodes (for the peripheral nervous system 29 ). 
Thus we shall consider a finite FHN system with N nodes and a Neumann boundary 
condition at the right end, un+i — Un- At the left end, we impose Uo(t) = 2 for 
< t < 0.05, and Uo(t) = for t > 0.05. The results corresponding to parameter val- 
ues d = 0.1 and a — 0.5 are depicted in Figures 3.1 (for which e = 0.006) and |6.2| (for 
which e = 0.003). The asymptotic theory predicts that fully developed FHN pulses 
(corresponding to N = oo) would have widths of I* s» 5 and I* w 10, respectively. 
The left boundary condition ensures that the membrane potential u n is excited during 
sufficient time, so that a wave is generated at the left end of the fiber. 

The excitation at the left boundary induces a wave front that propagates with 
velocity approximately given by C — c_ (0) /e along the finite fiber for the parameter 
values we consider. For example, C « 12.6 for e = 0.006, which is close to the 
numerically observed value of 10 in Fig. |6.1| . Similarly, C sa 25.22 for e = 0.003, which 
is close to the numerically observed value of 26 in Fig. S.2, If the fiber is long enough, 
a second wave front follows the first one and their mutual distance rapidly a ppro aches 
the asymptotic value I* (the number of nodes betw een fronts is 4 in Fig. 6T while 



the asymptotic theory predicts I* ~ 5; in Fig. 3.2, numerical observation confirms 



the asymptotic value I* « 10). The numerical solution of the finite FHN system 
shows that an eventually truncated FHN pulse comprising the two wave fronts and 
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the region between them is formed provided N is at least twice I*. Otherwise at best 



only the first wave front is shed at the boundary, as shown in Fig. 3.2(a). Pulses 
fail to propagate in fibers whose parameters fall in the propagation failure region, as 
discussed in Section |[ 

7. Conclusions. We have constructed stable pulses of the spatially discrete 
FitzHugh-Nagumo system by asymptotic methods. In a pulse, there are regions where 
the excitation variable varies smoothly separated by sharp fronts. These fronts are so- 
lutions of the discrete Nagumo equation with a constant value of the recovery variable. 
Their shape and speed can be calculated approximately near parameter values cor- 
responding to front propagation failure or near the continuum limit. For long times, 
their width is given by the only stable solution of a one-dimensional autonomous sys- 
tem. We have compared the asymptotic results with numerical solutions of the FHN 
system and analyzed different scenarios for failure of pulse propagation. Besides the 
classical scenario of scarce separation between time scales of excitation and recovery 
(large e as in the spatially continuous FHN system), propagation failure of fronts for 
the spatially discrete Nagumo equation provides a different mechanism of pulse propa- 
gation failure. Wave fronts and pulses can be generated at a boundary and propagate 
or fail to propagate along a finite FHN system. If the number of nodes is sufficiently 
large, the two wave fronts comprising a FHN pulse can be shed at the boundary and 
their separation rapidly reaches the value given by the asymptotic theory. This is so 
even if the fiber is too short to accommodate the slowly varying regions at the back 
of the second wave front of the pulse. In long enough fibers, a fully developed FNH 
pulse may be generated by an overthreshold stimulus applied during a short time at 
one end of the fiber. 
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